Astron. Nachr. / AN 999, No. 99, 999-[l004](2013) / DOI please setDOI! 



Dynamical stability of the Gliese 581 exoplanetary system 



Zs. Toth 1 * and I. Nagy 



2.3 



CO 

o 

CN 



CN 
Q-i! 

6 

CO 



CN 
> 
CN 
CN 
CO 

CN 
O 
CO 



13 



1 Department of Geosciences, University of Bremen, Klagenfurter StraGe, 23359 Bremen, Germany 

2 Department of Natural Science, National University of Public Service, Hungaria korut 9-11, 1 101 Budapest, Hungary 

3 Department of Astronomy, Eotvos Lorand University, Pazmany Peter setany 1/A, 1117 Budapest, Hungary 

Received 16 November 2012, accepted xx Month xxxx 
Published online xx 

Key words planetary systems - stars: individual (Gliese 581) - methods: numerical - methods: N-body simulations 

Using numerical methods we investigate the dynamical stability of the Gliese 581 exoplanetary system. The system is 
known to harbour four certain planets (b-e). The existence of another planet (g) in the liquid water habitable zone of the 
star is supported by the latest analysis of the radial velocity (RV) measurements. Vogt et al. (AN 333, 561-575, 2012) 
announced a 4- and a 5-planet model fitted to the RV curve with forced circular orbits. To characterize stability, we used 
the maximum eccentricity the planets reached over the time of the integrations and the Lyapunov Characteristic Indicator 
to identify chaotic motion. The integration of the 4-planet model shows that it is stable for even an inclination i = 5°, 
i. e. high planetary masses, on a longer timescale. We varied the orbital parameters of the innermost low-mass planet e, 
which quickly became unstable in earlier eccentric models. Planet e remained stable during the integration time, although 
the stable region around its initial semi-major axis and eccentricity is rather small. In the 4-planet model, we looked for 
stable regions for a fifth planetary body. We found extensive stable regions between the two super-Earth sized planets c 
and d, and beyond planet d. The Titius-Bode law and its revised version, Ragnarsson's formula applied to the Gliese 581 
planetary system both predict distances of additional planets in these stable regions. The planet Gliese 581 g would have 
a dynamically stable orbit, even for a wider range of orbital parameters. Since circular orbits in the models seem to be a 
too strong restriction and the true orbits might be elliptic, we investigated the stability of the planets as a function of their 
eccentricity, and derived dynamical constraints for the ellipticity of the orbits. 
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1 Introduction 



Based on high-precision ra dial velocity (RV) measurements 
acquired using the HAR PS dBonfils et all2005l: Mayor et ail 



cient to keep planet d from freezing out at the cold edge of 
the habitable zone (Sels is et al . 2007b- 



2009HUdrv et aD2007l) and HIRES ( Vogtetal 



20101) 



spec- 
trographs, an exoplanetary system was discovered around 
the red dwarf star Gliese 581 recently. The RV data sug- 
gest that the system harbours at least four, and maybe five 
or six planets: Gliese 581 e, b, c and d are certain, while 



the existence of f and g are disputed dForveille et al.ll2011 
Vogt etalJEoB) . 



The planetary system of Gliese 581 consists of a 
Neptune-mass planet on a 5.36-day orbit (planet b), two 
super-Earth sized planets (planet c and d) with orbital pe- 
riods of 12.9 and 67 days, and a low-mass planet closest to 
the star with an orbital pe riod of 3.15 days, all discovered in 



the HARP S data se t s bvlBonfils et alJ (120051) : iMay or et al 



d2009t) : lUdry et al.l d2007l) . respectively. The super-Earth 
planets in the system are located on the two edges of the 
liquid water habitable zone, and while the greenhouse ef- 
fect of the atmosphere would make planet c too warm and 
therefore unable to host liquid water, high concentrations of 
carbon dioxide or other greenhouse gases would be suffi- 



Combining the 4.3-yea r HARPS set with the 1 1 -year set 
of HIRES RVs. lVogt etall lfeoiO) announced the discovery 
of planet f and g orbiting the star with periods of 433 and 
36.6 days. Both planet f and g were indicated in the com- 
bined RV data sets using a Keplerian fit with (forced) cir- 
cular orbits. The quality of the fit could be improved only 
when the eccentricities on the 67-day and 36.6-day planets' 
orbits were allowed to float, with these two planets being 
in secular resonance. Gliese 581 g with a mass of ~ 3M ffi 
would be a rocky planet in the middle of the habitable zone 
of the system and might be habitable for a wider range of 
atmospheric conditions. 

Bayesian r e-analysis of bo th t he combined a nd individ- 
ual datasets bv lGregorvld201 lb and lTuomild201 II) confirmed 
only four clear planetary signals (planet e, b, c and d) and 
higher probability for the existence of planet f (with an or- 
bital period of ~ 400 days). Both studies found that the 
eccentricities for 3 of the 4 orbits are consistent with zero, 
the orbit of planet d is elli ptical (eg ^ 0.4), s imilarly to the 
earlier 4-planet solution of iMavor et al.l (12009b . 



Anglada-Escude et al.l (|2pi0) noted how solutions fit- 
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ting RV data sets with a single-planet eccentric orbit can 
hide two planets in circular resonant orbits. This is be- 
cause there is a degeneracy between the resonant and ec- 



® WILEY fj| ONLINE LIBRARY 



©2013 WILEY-VCH Verkg GmbH & Co. KGaA. Weinhcim 



1000 



Zs. Toth & I. Nagy: Dynamical stability of the Gliese 581 exoplanetary system 



centric solutions, as their Keplerian motion equations are 
identical up to the first order in the eccent r icity. I n a sub- 
sequent study, Anglada-Escude & Dawson] ( 20ld) showed 
that the first eccentricity harmonic of Gliese 581 d (~ 33.5 
days) coincides with a yearly alias of the newly reported 
planet g (~ 33.2 days), thus the high eccentricity of planet 
d can partially absorb the signal from planet g. Neverthe- 
less, based on statistical tests they concluded that the pres- 
ence of planet g is well suppo rted by the available RV data. 
Tadeu dos Santos et alj (120 121) similarly concluded that the 
existence of the 36-day planet g depends on the eccentric- 
ity of the 67-day planet d and its detection requires the as- 
sumption that all planets are on circular orbits. The signal 
of planet f was found, but only in the threshold of their con- 
fidence level with a period of ~ 455 days. 

released 



Forveille et al. 



d201 II) released an additional set of 
HARPS RV measurements and analysed the then total 7- 
year data set. Their 4-planet Keplerian-fits, with either 
fixed or freely floating eccentricities, revealed no signifi- 
cant residual signals, therefore no support for the two ad- 
ditional planets. The additional measurements revised the 
mass of planet d d own to 6 M m , making a rocky compo- 



sition more likely. IVogt et al.l (120 121) then re-analyzed this 
data set and warned that allowing the eccentricities to float, 
and in particular the eccentricity of planet e to rise, leads 
to instability and therefore to highly unphysical Keplerian 
models. None o f the N-body simulati ons of the eccentric 
Keplerian fit of Forveille et al. d201 lh remained dynami- 
cally stable on longer time scales due to high eccentricity of 
planet e. The new 4-planet all-circular interacting model of 
Vogt et al.l (120120 remains dynamically stable for 20 Myrs. 
Furthermore, it offers confirmative support for a fifth plan- 
etary signal near 32-33 days, which could be planet g at its 
36-day yearly alias period. 

Here we will present a numerical investigation into the 
dynamical stability of the Gli ese 581 system, s tarting from 
the 4 and 5-planet models of IVogt et al.l ( 2012h . The objec- 
tive of our study is to set dynamical constraints on the orbits 
of the latest RV fit, including planet g. Out of curiosity, we 
test how well the Titus-Bode law describes the Gliese 581 
planetary system and where it predicts additional planetary 
orbits. 



2 Titius-Bode law 

The Titius-Bode law (TBL) is an empirical formula de- 
scribing the planetary distances in the Solar System in the 
form of a power law. Historically, TBL described the plane- 
tary orbits from Mercury through Saturn, and correctly pre- 
dicted the orbits of Ceres and Uranus, but failed when Nep- 
tune and the Kuiper Belt objects were discovered. Although 
no simple physically relevant explanation of the TBL has 
been found so far, astrophysical s tudies point to physical 
background (see the summary of ICuntzll2012l) . and there 
has been attempts to investigate if exoplanetary systems fit 
Titius-Bode like laws, and to use them to predict yet undis- 



Table 1 Distances of the planets in the Gliese 581 plan- 
et ary system based on the Titius-Bode law and the formula 
of Ragnarsson ( 1995h . 



Planet 


a 

(AU) 




O.TBL 

(AU) 






ClRF 

(AU) 




e 


0.03 


n 


= — 00 


0.03 


m 


= -1 


0.022 


b 


0.04 


n 


= 


0.04 


m 


= 


0.04 


? 




n 


= 1 


0.05 








c 


0.07 


n 


= 2 


0.07 


m 


= 1 


0.074 


S 


0.13 


n 


= 3 


0.11 


m 


= 2 


0.148 


d 


0.22 


n 


= 4 


0.19 


m 


= 3 


0.221 


? 










m 


= 4 


0.294 


? 




n 


= 5 


0.35 


m 


= 5 


0.368 



covered exoplanets, e.g. in the 5-pl anet system of 55 Cancri 
dCuntzH20ia Poveda & Larall20 08). 

Based on the hypothesis of a geometrical planetary se- 
quence, we calculated the planetary distances in the Gliese 
581 system . Usin g the TBL and a revised formula by 
Ragn arsson (1 19951) . both adjusted to the known exoplanets 
of Gliese 581, we search for gaps for additional planets in 
the inner part of the system (> 0.4 AU). 

TBL states that the semi-major axis of the planets are 
given, in astronomical units, by the formula a n — 0.4 + 
0.3 ■ 2", n = -oo, 0, 1, 2, . . . Adjusted to the Gliese 581 
planetary system it can be written: 



a n = 0.03 + 0.01 • 2 n , n = -oo, 0, 1, 2, . 



(1) 



Ragnarsson's new empirical formula (RF) is more ac- 
curate in describing the Solar System, and is based on the 
symmetry in the logarithms of the semi-major axes around 
Jupiter. Assuming that the hot Neptune (planet b) in the 
Gliese 581 system plays the same dominating role as Jupiter 
in our Solar System, RF can be written: 

a m = a b (5/2) 2/3 |m| , to = -1, 0, 1, . . . (2) 



where ai, is the semi-major axis of Gliese 581b, param- 
eter m is the "jovicentric" planet number, which is negative 
for the inner planets (in this case for planet e). 

The distances of the planets calculated from Eq.(Q~|i and 
Eq.© are listed in Table Q] Neither of these formulas seem 
to describe perfectly the actual values, although both work 
for three out of the four planets relatively well, in fact it 
is remarkable how well RF fits for planet c and d. The 
TBL predicts three more possible planets in the system, but 
we question if at 0.05 AU, very close to the most massive 
planet, another stable orbit can exist. Between planet c and 
d, another planetary body is indicated at 0.11 and at 0.148 
AU by the TBL and the RF respectively. The next predicted 
orbit after planet d is at 0.29, ~ 0.36 AU. 

3 Dynamical models and methods 

We investigated the dynamical stability of the Gliese 581 
planetary system by using numerical integrations of the 
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Table 2 Astroce ntric, circu l ar, no n-interacting 4-planet 
orbital model of IVogt et al.l (120121) . P: orbital period, 



M m in — m sin i, where i is the orbital inclination, a: semi- 
major axis, e: eccentricity, /: longitude of the periastron. 



Planet 


p 




a 


e 


I 




(days) 


(M e ) 


(AU) 




(deg) 


e 


3.15 


1.84 


0.028 





138.5 


b 


5.37 


15.98 


0.04 





338.9 


c 


12.93 


5.4 


0.073 





175.2 


d 


66.71 


5.25 


0.22 





235.8 



planetary orbits. The model parameters w ere taken from the 
announced 4-planet and 5-planet fits of IVogt et al.l d2012l) 
and are presented in Table [2] and [3] for reference. For the 
mass of Gliese 581 we took 0.31 A/©. We made the assump- 
tion of coplanarity of all orbits, and in case of the 4-planet 
model we checked the stability for different orbital incli- 
nations (to the line of sight), otherwise and for the 5-planet 
model sin i = 1 (an edge-on system) was assumed. We con- 
fined our study area to the region of the four certain and 
the potential fifth planet, between 0.01 and 0.41 AU. In our 
study, we investigated the stability of the following models: 



- 5 -body problem consisting of the star and the 4 certain 
planets (b-e) to check the long-term stability of the 4- 
planet system, 

- restricted 6-body problem consisting of the star, 4 plan- 
ets (b-e) and a massless fifth planet to search for stable 
regions for an additional planetary body, 

- 6-body problem consisting of the star and 5 planets to 
verify the dynamical stability of the 5-planet syst em in- 
cluding planet g as presumed bv lVogt et al.l (120121) . 



The numerical integrations were performed using the 
method of Lie-integration with an adaptive step-size, which 
is a very fast and precise in tegration method due to the 
recurrence of th e Lie-terms (faanslmeie r & Dvorak 1984; 
Pal & Siilil 120071) . In order to characterize the stability of 
the models, we used the Lyapunov characteristic indicator 
(LCI), the finite time approximation of the maximal Lya- 
punov Exponent, and the maximum eccentricity method 
(MEM). The LCI is a well known chaos indicator of a dy- 
namical system, it estimates the exponential divergence rate 
of infinitesimally close trajectories in the phase space. The 
MEM provides information about the evolution of the or- 
bit during the integration time through the largest value of 
the eccentricity of the planets and indicates close encounters 
and escapes from the region of motion as well. 

4 Stability of the 4-planet system 

4.1 Basic stability of the 4-planet model 

The analysis of RV variations is able to constrain the mass 
(M) of exoplanets by a lower limit, since only the quantity 
M sin i is determined, where i is the unkown orbital incli- 
nation. Below a given value of i - or, equivalently above 



given masses for the planets - we expect the system to be- 
come unstable as the dynamical interactions can increase 
with lower inclinations, i. e. higher planetary masses. To 
check for which range of inclinations the system has this ba- 
sic stability, we integrated the 4-planet model (Tabled over 
1 Myr varying the inclination by Ai = 5° from i = 90° 
(edge on) to i = 5° (almost pole on). 

In all cases, the integrations show that the variations of 
the semi-major axes and eccentricities of the planets are not 
significant. The 4-planet system remains stable for 1 Myr 
down to i = 5° (almost pole-on) and co-planar orbits. Dy- 
namical stability with i > 5° suggests that the mass of each 
planet can be ~ 12.5 times of its minimum mass. For Gliese 
581 e, b, c, and d those limits are 21, 183, 62 and 60.2 Mm. 



In comparison, the 4-planet fit o flMayor et al. (2009) be- 
came unstable and planet e escaped from the system after a 
few Myrs already for i > 40°. Although their integration 
time was longer, the orbital parameters did not change sig- 
nificantly in our integration in case of smaller inclinations, 
which means considerably higher planetary masses and thus 
stron ger perturbations over the time of 1 Myr. The 4-planet 
fit of lVogt et alJ (12.012b in this sense can be considered more 
stable. 



4.2 Stability of planet e 

The low mass planet e has little effect on the stability of the 
more massive planets in the system, but an increase in its 
eccentricity quickly leads to orbital crossings with planet b 
or ejection from the system (as it w as already pointed out by 
Mayor et al.ll2Q09c IVogt et al.ll2012l) . To characterize planet 
e's orbit in the 4-planet model, we carried out integrations 
with the following grid of initial semi-major axis and eccen- 
tricity: 0.01 < a < 0.041 AU with steps Aa = 0.003 AU 
and < e < 1 with steps Ae = 0.01, for a time period of 
5000 P e , where P e is the orbital period of planet e. 



The maximum eccentricity that planet e reached over 
the time of the integration is plotted on the a — e stabil- 
ity map in Fig. [T] Planet e's orbit is stable (small e max , 
yellow color), although the vicinity around its initial a e = 
0.028, e e — values represents a stable island on the map, 
so relatively small perturbations could lead to instability and 
thus collision or ejection. Having varied only the eccentric- 
ity e e , the LCI became suddenly high above e e = 0.22 
(data not shown). This indicates chaos, and also points to 
the highest eccentricity for planet e's orbit to remain stable. 

We also performed integrations with different starting 
positions (I) for planet e varying I by 2° between I = 
— 360°. The LCI values remain low during the time of 
the integration for each case (data not shown), so planet e's 
stability does not depend on its orbital position in the model. 
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0.010 0.015 0.020 0.025 0.030 0.035 0.040 




Fig. 1 a-e stability map of planet e, its maximum eccen- 
tricity is plotted as a function of semi-major axis and eccen- 
tricity (with a starting value of 0). Yellow regions are stable, 
black regions denote escaping orbits. 



Fig. 3 LCI values of a fifth massless test planet for dif- 
ferent initial semi-major axis (a) and starting position (1). 
Chaotic orbits are marked with colors towards black, stable 
orbits towards yellow. 




Fig. 2 LCI values of a fifth massless test planet for differ- 
ent initial a-e. Stable regions are marked as yellow accord- 
ing to the LCI. 

5 Stability of the 5-planet system 

5.1 Stable regions for a fifth planet 

We integrated the 4-planet system with an additional mass- 
less body in order to find regions of possible orbital sta- 
bility for a fifth planet. The test planet's initial parame- 
ters were varied, in one set the semi-major axis between 
0.01 < a < 0.41 AU with steps of Aa = 0.001 AU and 
the eccentricity < e < 1 with steps of Ae = 0.01, and in 
another set for the same semi-major axis range and e = 0, 
the starting position < I < 360° with steps of Al = 2°. 
The time of the integrations was 10000 Ptest for calculating 
the LCI, where Pt es t is the orbital period of the test planet. 
The LCI values of each integration for the grid of different 
a, e are plotted on Fig. [2] for a, I on Fig. [3] The e max that 
the test planet reached over 5000 Ptest was calculated for 
the same grid of a, e, and for 8 different starting orbital po- 
sitions with steps Al = 45° between 1 = — 315°. The 
average of these 8 e rnax values are represented on one map 
on Fig.|U 

The orbits of a fifth fictitious planet are chaotic close to 
the innermost planets e and b, but extensive stable regions 
exist between the larger planets c and d, and outside planet 



0.4 



14.8 
4.0 

I 

2.4=1 

Icji 
o 
1 .6 I 

0.8 
0.0 



0.2 



0.3 



0.4 



Fig. 4 Maximum eccentricity of a fifth massless test 
planet. Each a, e point on the map represents an average of 
the maximum eccentricities that the test planet reached over 
the time of the integrations started from 8 different starting 
positions in the model. 



d, with relatively low eccentricities (Fig. |2). The starting 
position of the fictitious planet only plays a role in its sta- 
bility around the innermost planets and the orbit seems less 
chaotic in the outer part of the system (Fig. [3). The a — e sta- 
bility map with the maximum eccentricity essentially shows 
the same boundaries for stable regions, in fact the average 
e m ax values confirm that irrespective of the starting posi- 
tion of the test planet, its orbit will be stable in the region 
between the super-Earth planets and beyond them (Fig. 0J. 

There are small stable islands at the distance of planet d 
(a = 0.22 AU), where both the e max and the LCI have low 
values at the end of the integration (Fig. [3] and Fig.|4]l. This 
points to the possibility of Trojan planets orbiting the star in 
1 : 1 mean motion resonance with planet d. 

According to these stability maps, the additional plan- 
ets, predicted by the TBL and RF between planet c and d 
at 0.11 and 0.15 AU, and beyond planet d at 0.29 and 0.36 
AU, would have stable orbits with low eccentricities. 
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Table 3 Astroce ntric, circu l ar, no n-interacting 5 -planet 

d2012l) . P: orbital period, 
a: semi- 



orbital model of iVogt et al 
M m i n = m sin i, where i is the orbital inclination 



major axis, e: eccentricity, /: longitude of the periastron. 

Planet P M m i n a e I 
(days) (M e ) (AU) (deg) 



e 


3.15 


1.86 


0.028 





141.9 


b 


5.37 


16.0 


0.04 





338.4 


c 


12.93 


5.3 


0.073 





181.0 


a 


32.13 


2.24 


0.13 





55.3 


d 


66.67 


5.94 


0.22 





227.3 




Fig. 5 LCI values of the orbit of planet g as a function of 
eccentricity of the adjacent planets c and d. The yellow to 
red rectangular region indicates stable orbits, initial points 
in the dark regions result in chaotic motion. 



5.2 Stability of planet g 



The 5-planet model of lVogt et alj (12012) (as a reference see 
Table |3]l was integrated for 1 Myr to verify its long-term 
stability. All 5 planets remained stable during the long in- 
tegration time and their orbital parameters did not change 
significantly, so this 5-planet solution of the RV data can be 
considered dynamically stable. 

If planet g exists, it is a 2.2 M ffi planet orbiting be- 
tween two super-Earth sized planets in the middle of the 
star's classical liquid water habitable zone. To further char- 
acterize the stability of the orbit of this possibly rocky and 
habitable planet, we ran two more set of integrations. In one 
set, for a time of 5000 P g , the eccentricity of the adjacent 
planets c and d were varied until the point when their orbits 
cross each other. According to the LCI of planet g on Fig. 
|5j, planet g's orbit is not chaotic and thus can be considered 
stable in cases where the eccentricity of planet c stays below 
e c < 0.32 and for planet d e d < 0.28. 

Since the true orbital period and therefore semi-major 
axis of planet g is still uncertain, as well as planet d's real 
eccentricity could be greater than zero, in another set of 
integrations, also for a time of 5000 P g , we varied a g in 
the range between the neighbouring planets and until it 
would cross the orbit of planet c. Not surprisingly, as Fig. 
[6] shows, the more elliptic the orbit of planet d, the further 
away planet g has to be from it in order to remain stable. At 




Fig. 6 Stability of planet g as a function of planet d's ec- 
centricity and the semi-major axis. Yellow marks the stable, 
black end of the color bar the chaotic orbits according to the 
LCI. 



the largest distance, planet g's orbit would be still stable in 
case of ed ~ 0.55. 



6 The problem of orbital eccentricity 

In the process of fitting multi-planet models to precise RV 
data, it is a question of choice whether to allow the orbital 
eccentricities to float freely or to hold them fixed at zero. As 
the case of the Gliese 581 planetary system shows, freely 
floating eccentricities can lead to solutions that are unstable 
and therefore unphysical despite providing good fits to the 
RV data. On the other hand, forced circular orbits might be 
too strong a restriction. 

In the previous section we derived limits for the eccen- 
tricity of planets c and d with planet g in the system, now 
we will investigate the situation of the inner planets. We 
have assumed that given its higher mass planet b is likely to 
be more stable in comparison to the adjacent planets, there- 
fore we looked at the stability of planet e for a range of 
intial e e and e^, and the stability of planet c for a range of 
initial e& and e c (in this case planet e was left out of the 
system). The maximum eccentricity was calculated in both 
cases over a time period of 5000 P e and 5000 P c , respec- 
tively. The eccentricities ranged until the elliptical orbit of 
one planet crossed the other one's circular orbit. 

The maximum eccentricity on Fig. [7] indicates that 
planet e's stability is ensured when the adjacent planet b's 
eccentricity stays relatively low, e& < 0.075, and its own or- 
bital eccentricity below e e < 0.2, approximately the same 
upper limit as the one based on the LCI. Small stable re- 
gions exist with eccentricities above e e = 0.3 too, though 
these isolated configurations seem less likely to remain sta- 
ble for a longer time. Looking at the stability of planet c as 
a function of e\, and e c (Fig. [8]), it is apparent that planet c 
puts less constraint on the eccentricity of planet b (which is 
then mainly given by the smaller planet e). While in case of 
circular or near-circular orbits, the upper limit of planet b's 
eccentricity is e c < 0.32. 
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6.7 



1.3 



Fig. 7 Stability of planet e is characterized with its the 
maximum eccentricity as a function of e e and e&. Black de- 
notes unstable, whereas yellow stable configurations. 



Modelling the system in order to find a good fit to 
the RV curve was found to be quite complex, especially 
in regards to the question of how they treat the eccen- 
trities of the o rbits. The eccentric Keplerian model of 
Forveille et al. J 201 lb was very unstable, and while the 
models of Vogt et all ( 2012 ) with initial circular orbits are 
dynamically stable over a longer time period, the true or- 
bits might be elliptic as a consequence of dynamical inter- 
actions between the planets. Our simulations put dynamical 
constraints on the eccentricity of the planets, and these lim- 
its can in turn be used to help find the best stable model to 
fit the RV measurements. 

Acknowledgements. We greatly appreciate the help of William 
Brocas with the manuscript and we thank the anonymous reviewer 
for providing constructive comments. 




Fig. 8 Stability of planet c is characterized with its the 
maximum eccentricity as a function of e& and e c . Black de- 
notes unstable, whereas yellow stable configurations. 

7 Conclusions 

The long-term numerical integrations show that the latest 
4-planet all-circular model of the Gli ese 581 exoplanetary 
system presented bv lVogt et alj (120121) is dynamically stable 
for 1 Myr. The innermost planet e, which was found to be 
sensitive to small perturbations due to its low mass and the 
tight packing of the inner planets, remains stable, however 
the stable region on the a-e plane around it is rather small 
and its eccentricity cannot be larger than 0.2. 

A fifth planetary body in the 4-planet system could have 
a stable orbit between the two super-Earth planets c and 
d, and there is an extensive stable region for an additional 
planet beyond the orbit of planet d. The Titius-Bode law 
and Ragnarsson's formula applied to the Gliese 581 plane- 
tary system both predict planets at distances from the star, 
where large stable regions exist. 



The 5 -planet model of IVogt et alJ (120121) including 



planet g is also dynamically stable on a longer timescale. 
The existence of planet g is fully supported from the dy- 
namical point of view. In fact a low-mass planet like planet 
g can be considered stable for a wider range of orbital con- 
figurations in the region between the two super-Earth plan- 
ets. 



References 

Anglada-Escude, G., Lopez-Morales, M., Chambers, J. E.: 2010, 

ApJ709, 168-178 
Anglada-Escude, G., Dawson, R. I.: 2010 ApJL, submitted; 

larXiv:1011.0T86V 2 [astro-ph.EP] 
Bonfils, X., Forveille, T., Delfosse, X., et al.: 2005, A&A 443, 15- 

18 

Cuntz, M.: 2012, PAS J 64, 73 

Hanslmeier, A., Dvorak, R.: 1984, A&A 132, 203-207 

Forveille, T, Bonfils, X., Delfosse, X., et al.: 201 1, A&A, submit- 
ted; [^^11092505}' 1 [astro-ph.EP] 

Gregory, P. C.: 2011, MNRAS 415, 2523-2545 

Mayor, M., Bonfils, X., Forveille, T., et al: 2009, A&A 507, 487- 
494 

Pal, A., Siili, A.: 2007, MNRAS 381, 1515-1526 

Poveda, A., Lara, P.: 2008, RMxAA 44, 243-246 

Ragnarsson, S-I.: 1995, A&A 301, 609-612 

Selsis, F, Kasting, J. F, Levrard, B., Paillet, J., Ribas, I., Delfosse, 

X.: 2007, A&A 476, 1373-1387 
Tadeu dos Santos, M., Silva, G. G., Ferraz-Mello, S., 

Mtchtchenko, T. A.: 2012, CeMDA 113, 49-62 
Tuomi, M.: 2011, A&A 528, L5 

Udry, S., Bonfils, X., Delfosse, X., et al.: 2007, A&A 469, 43-47 
Vogt, S. S., Butler, R. P., Rivera, E. J., Haghighipour, N., Henry, 

G. W., and Williamson, M. H.: 2010, ApJ 723, 954-965 
Vogt, S. S., Butler, R. P., Haghighipour, N.: 2012, AN 333, 561- 

575 



©2013 WILEY- VCH Verlag GmbH & Co. KGaA, Wcinhcim 



www.an-journal.org 



